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Abstract Extended Chebyshev spaces that also comprise the constants represent large families of functions that 
can be used in real-life modeling or engineering applications that also involve important (e.g. transcendental) 
integral or rational curves and surfaces. Concerning computer aided geometric design, the unique normalized 
B-bases of such vector spaces ensure optimal shape preserving properties, important evaluation or subdivision 
algorithms and useful shape parameters. Therefore, we propose global explicit formulas for the entries of 
those transformation matrices that map these normalized B-bases to the traditional (or ordinary) bases of the 
underlying vector spaces. Then, we also describe general and ready to use control point configurations for the 
exact representation of those traditional integral parametric curves and (hybrid) surfaces that are specified 
by coordinate functions given as (products of separable) linear combinations of ordinary basis functions. The 
obtained results are also extended to the control point and weight based exact description of the rational 
counterpart of these integral parametric curves and surfaces. The universal applicability of our methods is 
presented through polynomial, trigonometric, hyperbolic or mixed extended Chebyshev vector spaces. 
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1 Introduction 

Normalized B-bases (a comprehensive study of which can be found in [24] and references therein) are normalized 
totally positive bases that imply optimal shape preserving properties for the representation of curves described 
as convex combinations of control points and basis functions. Similarly to the classical Bernstein polynomials of 
degree n £ N - that in fact form the normalized B-basis of the vector space of polynomials of degree at most n 
on the interval [0,1], cf. [9] - normalized B-bases provide shape preserving properties like closure for the affine 
transformations of the control points, convex hull, variation diminishing (which also implies the preservation of 
convexity of plane control polygons), endpoint interpolation, monotonicity preserving, hodograph and length 
diminishing, and a recursive corner cutting algorithm (also called B-algorithm) that is the analogue of the 
de Casteljau algorithm of Bezier curves. Among all normalized totally positive bases of a given vector space 
of functions the normalized B-basis is the least variation diminishing and the shape of the generated curve 
more mimics that of its control polygon. Important curve design algorithms like evaluation, subdivision, degree 
elevation or knot insertion are in fact corner cutting algorithms that can be treated in a unified way by means 
of B-algorithms induced by B-bases. 

Curve and surface modeling tools based on non-polynomial normalized B-bases also ensure further advan¬ 
tages like: possible shape or design parameters; singularity free exact parametrization (e.g. parametrization of 
conic sections may correspond to natural arc-length parametrization); higher or even infinite order of precision 
concerning (partial) derivatives; ordinary (i.e., traditionally parametrized) integral curves and surfaces can ex¬ 
actly be described by means of control points without any additional weights (the calculation of which, apart of 
some simple cases, is cumbersome for the designer); important transcendental curves and surfaces which are of 
interest in real-life applications can also be exactly represented (the standard rational Bezier or NURBS models 
cannot encompass these geometric objects). Moreover, concerning condition numbers and stability, a normal¬ 
ized B-basis is the unique normalized totally positive basis that is optimally stable among all non-negative 
bases of a given vector space of functions, cf. [24, Corollary 3.4, p. 89]. These advantageous properties make 
normalized B-bases ideal blending function system candidates for curve (and surface) modeling. 
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Besides their interest in the classical contexts of CAGD and approximation theory, normalized B-bases and 
their spline counterparts have also been used in isogeometric analysis recently (consider e.g. [21] and references 
therein). Compared with classical final element methods, isogeometric analysis provides several advantages 
when one describes the geometry by generalized B-splines and invokes an isoparametric approach in order 
to approximate the unknown solutions of differential equations (e.g. of Poisson type problems) or Dirichlet 
boundary conditions by the same type of functions. 

Let n > 1 be a fixed integer and consider the extended Chebyshev (EC) system 

= {V>n,i («) : U £ [a,/9]}" =0 , ipn.O = 1, —00<01<I3<00 (1) 

of basis functions in C n ([a, /?]), i.e., by definition [15], for any integer 0 < r < n, any strictly increasing 
sequence of knot values a < mo < ui < • • • < u r < (3, any positive integers (or multiplicities) such that 

Ylk=o m k = n + 1, and any real numbers k there always exists a unique function 

n 

f-=Yl := := span , A n>i e R, i m 0,1,..., n (2) 

i =0 

that satisfies the conditions of the Hermite interpolation problem 

f W (ttfc) = £k,l, l = 0,1, ■ ■ •, rn k - 1, k == 0,1,..., r. (3) 

In what follows, we assume that the sign-regular determinant of the coefficient matrix of the linear system 
(3) of equations is strictly positive for any permissible parameter settings introduced above. Under these 
circumstances, the vector space of functions is called an EC space of dimension n + 1. In terms of zeros, 
this definition means that any non-zero element of §n’^ vanishes at most n times in the interval [a,/3]. Such 
spaces and their corresponding spline counterparts have been widely studied, consider e.g. articles [17,10,22, 
23,18,16,6,19,13,8,20] and many other references therein. 

Hereafter we will also refer to as the ordinary basis of §“ ,/3 . Using [11, Theorem 5.1] it follows that the 
vector space also has a strictly totally positive basis, i.e., a basis such that all minors of all its collocation 
matrices are strictly positive. Since the constant function 1 = <p ni o £ §n ,/3 , the aforementioned strictly positive 
basis is normalizable, therefore the vector space also has a unique non-negative normalized B-basis 

Bn jj = {bn,i (u) :u£ [a, /3]}? =0 (4) 

that besides the identity 

n 

^2b n ,i (u) = 1, Vi* € [a,P] (5) 

2 = 0 

also fulfills the properties 

bn ,0 (o) = bn,n (/3) — 1, 
b n'} («) = 0, .7 = 0,...,*— 1, b^\ (a) > 0, 

b n\ 09) » o, 3 = 0,1,... , n - 1 - i, (-1)"-' (13) > 0 

conform [11, Theorem 5.1] and [22, Equation (3.6)]. 

Using the normalized B-basis of §“’^, one of our objectives is to provide explicit closed formulas for 
the control point based exact description of integral curves that are specified with coordinate functions given 
in traditional parametric form in the ordinary basis of the same vector space. Based on homogeneous 

coordinates and central projection, we also propose an algorithm for the control point (and weight) based exact 
description of the rational counterpart of these ordinary integral curves. Results will also be extended to the 
exact representation of families of (hybrid) integral and rational surfaces that are exclusively given in each of 
their variables by using ordinary EC basis functions of the type (1). 

To the best of the author’s knowledge, the coefficient based exact representation of ordinary (rational) 
functions, curves and surfaces by means of the (rational or spline counterpart) of the normalized B-basis of an 
arbitrary EC space (that also comprises the constant functions) was not considered in such a general unified 
context. Without providing an exhaustive survey, so far the presented problem appears in the literature for 
example in case of conversion algorithms related to Bernstein polynomials, monomials and the classical families 
of orthogonal Jacobi, Gegenbauer, Legendre, Chebyshev, Laguerre and Hermite polynomials [4,2] in special 
lower dimensional vector spaces (e.g. in [35,18,5,7,27]); in case of conical and helical arcs, of catenaries, of 
patches on all types of quadrics and of helicoidal surfaces (e.g. in [25,16]); of certain (rational) trigonometric 
curves of arbitrarily finite order like epi- and hypotrochoidal arcs [30] , or segments of offset-rational sinusoidal 
spirals, arachnidas and epi spirals [31]; or more recently, in case of arbitrary trigonometric and hyperbolic 
(rational) polynomials, curves, (hybrid) surfaces and volumes of finite order [28]. 

The rest of the paper is organized as follows. Section 2 lists our main results, namely it describes closed 
formulas for the basis transformation that maps the normalized B-basis B“’' 3 of the vector space §“ ,/3 to 
its ordinary basis J~n'' j and also specifies control point configurations for the exact representation of certain 
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large classes of integral and rational curves and surfaces that are specified in traditional parametric form by 
means of ordinary bases like Although the presented results are mainly of theoretical interest, Section 

2 also studies the computational complexity of the proposed basis conversion formulas and - compared with 
alternative cubic time numerical methods like curve interpolation or least square approximation - points out 
that these can more efficiently be implemented up to n = 15. Section 3 emphasizes the universal applicability 
of the general basis transformation described in Section 2 with examples that can be compared to presumably 
already existing results in the literature. This section considers EC vector spaces of functions that may be 
important in computer aided geometric design, in engineering, in (projective) geometry, in (numerical) analysis 
or in approximation theory. The proofs of all theoretical results stated in Sections 2-3 can be found in Section 
4. In the end, Section 5 closes the paper with our final remarks. Based on the general context of the manuscript, 
Appendix A recalls the classic transformation matrix that maps the Bernstein polynomials of degree n to the 
corresponding ordinary power basis of the vector space of traditional polynomials, while Appendix B provides 
implementation details by means of a simple Matlab example. 


2 Main results and remarks 


At first, we provide explicit formulas for the transformation of the normalized B-basis of the vector space 
to its ordinary basis 

Theorem 2.1 (General basis transformation) The matrix form of the linear transformation that maps the nor¬ 
malized B-basis Bnto the ordinary basis is 


[Pn,i ( u ) ] i=0 ~ [**j]i=o, j=0 ' i bn d ( M ) ]»=0 1 € ^ ’ 

where tg j = 1, j = 0,1,..., n and t™ 0 = i (a), t" n = iprM (0) , * = 0,1,..., n, while 


^i,j — ’Pn^i (°0 


h nl "‘ ! .=1 b n) («) 

j—r—1 

+ E (-tf E 

1= 1 r<k\ <k- 2 <...<k f <j 

*-*•*. n ' i=h2 .L1J • 


bn}r («) + 


(a) b™ (a) b™ (a)... (a) (a) 


u(k 2 


h( k s) 


(k t ) 


dj) 


(9) 


( 10 ) 


+ 


Vn?i ( a ) 

b nJ («) ’ 


j-1 , n ( r ) 


ti,n—j — P’n^i (/^) 


E 


,(r) 


b n,l-j (0) ^1 b n,n-r (0) 


b U) (0) + 

n : n—r \rJ 1 


(ii) 


CU ( 0 ) b n,n- kl ( 0 ) b n,n-k 2 ( 0 ) ■ ■ ■ CU-x ^ b Z~kt W 


b n'n-k, ( 0 ) b n,n-k 2 (?)■■■ b n,n-k e ( 0 ) 


j—r—1 

+ E (-i ) 4 E 

1=1 r<k 1 <k 2 <...<ke<j 

Vn} OT . , „ . , „ » | 

+ ,= ^. n ' , = 1 ' 2 .LaJ- 

n,n—j \' 0 


u(k 3 ) 


(kt) 


k(3) 


u(k 2 ) 


(kt) 


Remark 2.1 (Evaluation) If in formulas (10) or (11), for some I = 1, 2,..., j — r — 1 (with r = 1, 2,..., j — 1 and 
j = 1,2,..., \jl\) there exist no integers ki, & 2 ,... ,kg such that r < k\ < &2 < ■ • ■ < kg < j then, by convention, the 
summation corresponding to I equals 0. If n = 2z > 2, then for j = z one can evaluate the entries [tj, z ]” =1 of the 
middle column by using either of these formulas, since the zth coefficients of the ordinary basis functions (1) in the 
normalized B-basis (4) are unique. 


Except some special but important cases, in general, one does not know the closed form of the normalized 
B-basis (4) of §“’^. In case of EC spaces of traditional, trigonometric or hyperbolic polynomials of finite degree 
we have explicit closed formulas cf. [9], [29] and [32], respectively; in case of a special class of mixed (e.g. 
algebraic trigonometric, algebraic hyperbolic, or both trigonometric and hyperbolic) EC spaces these functions 
appear in recursive integral form cf. [20] and references therein; while the most general (determinant based) 
formulas that can be applied in such spaces was published in [22]. Thus, concerning the evaluation of (10) and 
(11), in general, one can differentiate the formulas presented in [22, Theorem 3.4, p. 658] in order to calculate 
the higher order derivatives of the normalized B-basis functions (4) at the endpoints of the interval [a,0]. 
Namely, by using the function 


0(«) ;= [<P„,1 (u) (fin,2 (u) 1 fin,n(u)] , U € [a, 0\ , 

one has to substitute the parameter values u = a and u = 0 into the derivative formulas 
det [0«(/?) </> {n - 1) (0) cf ij) (u)] 


b n,0 ( u ) = 


det[^(0) ... ^-^(0) 0(a)-0(/J)] 


( 12 ) 
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U) , , = detf^Ca) ■■■ 0("~ 1 >(a) 0^ (u) ] 

n ’ n(U> det[0«(a) ••• 0 ( 0 ) - 0 (a) ] ’ 

6 a) = det[0( 1 )( Q ) ■■■ ^-^(q) 0« (a) 0^ (/?) ■■■ 0C"-*" 1 ) (0) 0^> (/?) ] 

“ det[0W(a) 0 (i_1) (a) 0(0)-0(a) 0« (0) {P) 0 (n ” l) (0) ] 

det[0(0 )-0(a) 0^(a) ••• 0^ (a) 0W> (u) 0 (1) (/3) (0) ] 

' det [ 0(0)-0(a) 0 (1) («) 0 (i_1) («) 0 W («) 0 (1) (jS) ■■■ 0 (n " i_1) (/3)] 


for all i = 1,2,..., n — 1 and j = 1,2,..., n. However, as it is also mentioned in [22], these general relations 
are difficult and computationally expensive to evaluate even in the most simple cases for either arbitrarily 
big or general values of the order n. Therefore, Section 3 provides explicit closed formulas for the required 
endpoint derivatives in several special cases. Due to properties (7) and (8), these expressions should only be 
used whenever one does not know the exact value of the required endpoint derivatives. 

Another core result of the current section is presented in the next statement which is an immediate corollary 
of Theorem 2.1. 


Corollary 2.1 (Exact description of ordinary integral functions) Let {Ai}”_ 0 be real numbers and consider 
the linear combination 

n 

C ( U ) = 52 X iTn,i ( U ) ’ U£ I®’ £] ( 15 ) 

i =0 

of ordinary basis functions. Then, we have the equality c{u) = J2'j=oPjb n ,j ( u ) , V« £ [a, P\, where Pj = x P%,j, 

j = 0, l,...,n. 

Based on Corollary 2.1, the exact description of ordinary integral curves as convex combinations of control 
points and normalized B-basis functions (4) is presented in the next theorem. 

Theorem 2.2 (Exact description of ordinary integral curves) The ordinary integral parametric curve 

n g 

c(u) = (u) , u e P] > = [ A i] _ eR 6 , 5> 2 (16) 

2=0 _ 

of order n can be written as an EC B- curve 

n ^ 

c ( u ) = 52 P J b n,j ( u ), G [a, P\ , Pj = [pj] G It 5 , (17) 

3=0 

of the same order, where pj = X^”=o X i^i,j, j~ 0,l,...,n, £ = 1, 2,..., 5. 

Using tensor products of convex combinations of type (17), one can also exactly describe large families of 
surfaces as it is specified in the following theorem. 

Theorem 2.3 (Exact description of ordinary integral surfaces) Let 

= {<Pn r ,i r (u r ) ■ Ur G [ot r , Pr]}^ =Q , Tn r ,0 = 1, V — 1,2 

be two ordinary EC bases of some vector spaces §“0' Sr ' of functions and also consider their unique normalized 13- 
bases Bnp^ r = {bn r ,j r ( u r ) : u r G [a r , pr] }J r _ 0 , r* — 1,2. Denote by [tif j r ]^=Q r j r=0 the regular square matrix that 
transforms Bn(’^ r to and consider the ordinary integral surface 

s(u) = [s 1 (u) s 2 (u) s 3 (u)] T G R 3 , u = [u r ] 2 r= 1 G [ai,0i] x [ 02 ,02] (18) 

of order n = [n r ] r=1 , where 


(j£ 2 / n. 


S' ( u ) = 52 II ( 52 X i’rVn r ,ir ( u r ) ] , > 1, £=1,2,3. 

(=lr=l \i r =0 J 

Then, the surface (18) can be written in the tensor product form with the EC B-surface 

ni n 2 

s ( u ) = 52 52 P h,h b nuh ( u i) b n 2 ,j 2 (u 2 ), Vu = [u r ] 2 =1 G [ai,0i] X [a 2 ,0 2 ] 

ii=o j 2 =o 

of the same order, where the vectors Pj 1 ,j 2 = [pj L j 2 ]g=i G R 3 form the control net defined by coordinates 

<J£ 2 


(19) 


( 20 ) 


V U J ,C ._ U.C, 
C—1 r —1 i T =0 


pi, j2 = 52 n pf ’ Pjf := 52 A t C C,A . 3r = o, 1,..., nr, r = 1, 2, £ = 1, 2,..., a t , £ = 1, 2, 3. 


( 21 ) 
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Remark 2.2 (Exact description of ordinary integral volumes) Naturally, Theorem 2.3 can easily be extended 
to the control point based exact description of those tri- or higher variate integral multivariate surfaces (volumes) that 
are specified in traditional parametric form with coordinate functions described as sums of separable products of linear 
combinations of the type (15). 

If the denominator of the rational counterpart of the ordinary integral curve (16) is strictly positive, then, 
by means of control points and non-negative weights of rank 1, one can also exactly describe ordinary rational 
curves as it is illustrated in the steps of the next algorithm. 

Algorithm 2.1 (Exact description of ordinary rational curves) Consider in R 5 the rational curve 

C (W) = c *+i ( u ) V M t=i ’ U G [a ’ ^ (22) 

given in ordinary parametric form, where 

n 

c (“) = 53 X *Tn,i («), £=1,2,...,5 + 1, c S+1 (u) >0, Vwe [a, 0]. 

2 = 0 


Using the rational counterpart of EC B-curves (17), the process that provides the control point and weight based exact 
representation 

o w jPjbn,j (w) 


:(u) = 


Er=0 w rbn,r (u) 


Vu € [a, 0\ 


consists of the following steps: 

— apply Theorem 2.2 to the higher dimensional pre-image c p ( u ) = | '</ ( u ) 


points pj = 


jl 


5+1 


5+1 
£= 1 ' 


(23) 


u G [a,/3], i.e., compute control 




, j = 0,1,..., n for the exact description of c p in the pre-image space 


— project the obtained control points from the origin 0^_|_i £ R' 5 ^ 1 onto the hyperplane a:^ 1 = 1 that results in the 


control points p , = 3 [p(l 

J „5 + l L J \t=l 


^ and weights Wj = pj +1 needed for the rational representation (23); 


— the above generation process does not necessarily ensure the non-negativity of all weights, since the last coordinate 
of some control points p p in the pre-image space R' 5-1 " 1 can be negative; if this is the case, one should elevate 

the dimension (and consequently the order n of the normalized B-basis Bn'^) of the underlying EC space with an 
algorithm that generates a sequence of control polygons in R^^ 1 that converges to c p which, by definition, is a 
geometric object of one branch that does not intersect the vanishing plane x s+1 = 0, since the (5 +- l)th coordinate 
of all its points are strictly positive; therefore, by using proper dimension elevation methods, it is guaranteed that 
exists a finite and minimal order for which all weights are non-negative. 


Remark 2.3 (About the last step of Algorithm 2.1) If the pre-image c p of (22) is described as B-curves of 
type (17) by means of the normalized B-bases of the EC spaces C 8“^, then 


n n+1 

C P («) = Pj b n,j ( u ) = 53 Pl,j b n+l,j (u), Vu G [a, 0] , 

3 =0 3=0 


where pf 0 = p$, pf n+1 = p£, while pfj = (1 - + £jpJ for some real numbers £ (0,1), j = 0,1,... , n. 

Iterating this comer cutting based representation of c p in the normalized B-bases of the nested EC spaces §n ,/3 C 
§^3+ C ... C C one obtains a sequence of control polygons which converges to a Lipschitz-continuous 

limit curve [3] that, in general, does not necessarily coincide with c p . As it is pointed out in a unified manner in [28, 
Remark 2.3, p. 76], in case of vector spaces of finite order trigonometric/hyperbolic polynomials, the sequence of order 
elevated control polygons always converges to the curve generated by the first term of the sequence. In case of traditional 
polynomials of finite degree, one can use the well-known degree elevation techniques of (rational) Bezier curves. 
However, in general , the initial ordinary basis Tn'^ can iteratively be appended by new linearly independent functions 
in infinitely many ways and not every choice of functions leads to a sequence of order elevated control polygons that 
fulfills the desired convergence property, e.g. in EC Muntz spaces a recent characterization of the required convergence 
of the dimension/order elevation process can be found in [1], In order to illustrate the last step of Algorithm 2.1, Fig. 
1 shows two different control point configurations for the exact representation of the rational trigonometric curve 


c(u) 


1 

I - I sin (2u) 


77 g (sin (u) + cos (u)) 
| cos(2u) 


u e 



(24) 


by means of second and fourth order normalized trigonometric basis functions of type (29). 
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order trigonometric normalized B-basis functions of type (29) do not ensure the non-negativity of all weights 
required by the rational representation (23), since the third coordinate W 2 of the control point pf? is strictly 
negative. Such pathological cases usually violate the convex hull property of (23). Increasing the order of the 
applied normalized B-basis functions from 2 to 4, all control points needed for the exact representation of the 
pre-image c p will be above the vanishing plane, i.e., in this case there exists a minimal order starting from 
which the applied dimension elevation technique ensures the non-negativity of all weights. 


Remark 2.4 (Exact description of ordinary rational surfaces) The steps of Algorithm 2.1 can easily be ex¬ 
tended to the control point based exact description of those ordinary rational surfaces 


; ( u ) = TH ( U ) ( U ) s3 ( U )] T G R3 > U = Klr= 1 e [ai,Pl\ X [<X2,02] , 

s 


(25) 


in case of which 


<Jt 2 


! W = ^n E A lfg>n T ,ir K) ), ae> 1, € = 1, 2,3,4 

(=lr=l \v=0 


and s 4 (u) >0, Vu € [ai,/?i] x [ 0 : 2 , # 2 ]- 

Theorem 2.4 (Computational complexity) Provided that endpoint derivatives ( a ) i Tn \ iP) > ( a ) ; 


h LfJ 


^n \ 2 j=o are ca ^ cu ^ a l e d an d stored in advance in permanent lookup tables for a fixed value of n > 1, the 

total number of floating point operations (or shortly flops) that have to be performed for the evaluation of formulas 
(10) and (11) is 


k ( n) 


2 UJ+ 1 (L^J -3) +2 L^J +6 + 2n [§J ([fj + l) , n (mod2) = 1, 

2?- 1 _ 10 ) +n + 5+ rA n (mod2) = 0. 


(26) 


Using the normalized B-basis the control point based exact description of the ordinary integral curve 

(16) can also be imagined either as a curve interpolation problem or as the least square approximation of 
the considered curve. Both of these alternative numeric methods can be reduced to the solutions of 5 systems 
of linear equations that determine the unknown coordinates of the control points [Pj] n =0 appearing in the 
EC B-curve representation (17). Such methods also depend heavily either on the choice of the interpolation 
conditions or on the applied quadrature formulas used for the approximation of those integrals that appear in 
the equivalent quadratic form of the L 2 distance function used in case of least square approximations — but 
let us neglect both the floating-point round-off errors and the computational cost (i.e., the number of flops) of 
the regular main matrices of the size (n + 1) x (n + 1) of these alternative methods, and let us compare the 
exponential computational cost (26) to the total work of a possible algorithm that efficiently solves J systems 
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of linear equations. The fastest currently known matrix inversion algorithms [12] and [34] are based on the fast 
matrix multiplication algorithm of Strassen [33] and have an asymptotic cost of order O n + i) 2 - 376 ^ and 

O ({n + l) 2 ' 373 j. respectively, instead of O [{n + l) 3 ^ of traditional matrix inversion algorithms based e.g. on 
LU decomposition. However, if one estimates how large n+1 has to be before the difference between exponents 
2.373|2.376 and 3 is substantial enough to outweigh the bookkeeping overhead, arising from the complicated 
nature of the recursive Strassen algorithm, one finds that LU decomposition is in no immediate danger of 
becoming obsolete [26, p. 108]. The fast matrix multiplication/inversion algorithm of Strassen is typically used 
for n + 1 > 99, thus it is not practical to implement it for modeling purposes, since it is very unlikely that one 

would describe an arc of an ordinary integral curve with more than 99 control points. 

In case of a regular main matrix of the size (n + 1) x (n T 1), the number of flops performed by a numerical 
curve interpolation or least square approximation method based on LU decomposition is 

klu (n,5) = | ( n + l) 3 ~ \ ( n + 1 ) 2 - ^ (n + i) + (2 (n + l) 2 - (n + 1)) S (27) 

which covers the cost of computing of all multipliers, of all row operations and of all forward and backward 

substitutions as well. Naturally, as n tends to infinity, the growth rate of the exponential cost function (26) is 
substantially bigger than that of the cubic one (27), however if n £ {1, 2,... , 15} then (26) is less than (27) 
as it is illustrated in Fig. 2, i.e., compared with other cubic time numerical algorithms, the proposed general 
basis transformation can more efficiently be implemented up to 16-dimensional EC spaces despite the seemingly 
complicated nature of formulas (10) and (11). Considering that, in practice, curves and surfaces are mostly 
composed of continuously joined lower order arcs and patches, even a sequential but clever implementation of 
Theorem 2.1 can be useful in case of real-life applications. Nevertheless, if n > 15 then the presented results 
are mainly of theoretical interest. 



Fig. 2: Logarithmic scales of computational costs (26) and (27) for different values of n > 1 and 8 > 1. 


3 Examples 

This section applies closed formulas (10) and (11) in case of different vector spaces of functions that can be 
spanned by ordinary EC bases of the type (1). Our intention is only to emphasize the global applicability of the 
general basis transformation described in Theorem 2.1 with examples that can be compared to possible already 
existing results in the literature. Formulas (10) and (11) depend on the higher order endpoint derivatives of 
the ordinary and normalized B-basis of the underlying vector space. The following subsections specify these 
values in case of vector spaces of functions that may be important in many areas of applied or computational 
mathematics. We consider several reflection invariant EC spaces, since in practice usually one uses unbiased 
or symmetric systems of basis functions that also provide some computational advantages. Naturally, general 
formulas (10)-(11) are valid in not necessarily reflection invariant EC spaces as well. 


3.1 Trigonometric polynomials 


Let a = 0 and /3 £ (0,7r) be fixed parameters and consider the ordinary basis 

856 {^ 2 n,o (w) = 1 , {<p 2 n, 2 i-i (u) = sin (iu) , <p 2n , 2 i (u) = cos (iu)}" =1 : U £ [0, P\} (28) 

of trigonometric polynomials of order at most n (degree 2 n). Using the results of [29], the normalized B-basis 
of the vector space = \ J~ 2 n / can linearly be reparametrized into the form 
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where 


0 /3 

C 2n,i C 2n,2n—i 


LiJ 


sin 2 ™ ( f ) r =o 


E 


" '2coefg 

i — r I \ r I \ \ 2 


i = 0,1,..., n 


(30) 


are symmetric normalizing coefficients. It is obvious that 

v4n,2i-l (°) = iJ V2n,2i -1 , ¥>2n,2i-l (Z 3 ) = *V2n,2i-l ^ + y J ’ 3 = 0, i = 1,2...... 

<P2n,2i (°) = i3 V2n,2i ^y) > ^2n,2i (@) = $ V2n,2 i (p + y^j , j - 0,1,.. .,71, i = 1,2, ...,71, 

&2n,0 (0) = l*2n,2n (/?) = 1, &2n,i (0) = b 2n ,2n-i (P) =0, 7 = 1,2, ... , 2n, 

while the higher order derivatives {ly/l t (0), 4nt (/ 3 )}i=o”j=i are specified by the next theorem. 

Theorem 3.1 (Trigonometric endpoint derivatives) For arbitrary derivative order j = 1,2,... ,n we have that 

(0) 


^2n,2r+l 1 


r P 

L 2n,2r+1 


71 / — T*—1 T* / \ / \ 

^ E E(-i)" +1 -‘-'f 2(n_ 7 1) + 1 )f 2 T)' 

k=o e=o V / V / 

((n — k — ty — (71 — k — 2r + £ — l) 3 j cos ^(2 (n — r — k) — 1) ^ — y ^ 


(31) 


/or all r = 0,1,..., n — 1 and 


b (J> (0) (2n-2 r) r— 1 

u 2n,2r ) _ \ n-r > / 1 \i 

—J “ 2 2 "' 1 ^ 1 


¥=0 


—l / 2r \ . . ,• / jn 

1 1 (r — ly cos ( 


(32) 


+ ^k^T E (- 1 )" r k [ 2 ^ n k ^ ) (n-r ~k) J cos ((n-r-k)f3- J7r 


k =0 
n—r—1r —1 


+ 92^=1 E E(-^ 


2 2n-l \ fc / l £ 

fc=0 £=0 


-fc-< (2(n-r)\ [2r 



J7T 

2 


• ((n — k — (.y + (n — k — 2r + l) p ) cos ^(n — r — k)/3 - 

for all r — 0,1,..., n. At the same time b ^ i (ft) = (—l) 3 b^n 2 n-j (0) i * = 0,1,..., 2n, j = 0,1,..., n. 

Example 3.1 (Second order trigonometric polynomials) Consider the ordinary basis 

= {¥>4,0 (u) = 1, ¥>4,1 (u) = sin (u) , ¥>4,2 (u) = cos (u) , 

¥> 4,3 (u) = sin (2 u), ip 4,4 (u) = cos (2u) : u £ [0, /?]} , /3 G (0,7r) 

o/ £/ie vector space of trigonometric polynomials of order at most two (or degree 4) and its normalized B-basis 


>3°' 13 = \ 6 4 ,i («) = cf i sin 4 1 ( ) sin * (|) : w G [0, 


(33) 


where c P - r P - 1 c P - c P - 4coa (f) r p - 2 + 4cos2 (l) 

where c 4 0 - c 4 4 - —, c 41 - c 4 3 - —c 4:2 - —^ , 


and 


= 1, 64.1 

(0) 

= 0, 




1*4,2 

(0) 

= 0, 


1*4,3 (0) 

= 0, 



1*4,4 (0) 

= 0, 64,1 

(0) 

= 0, 




1*4,2 

(/l) 

= 0, 


1*4,3 iP) 

= 0, 



b 4 ,4 (/3) 

6 (1) 

4,1 

(0) 

J 

< 

h sir 



6 (1) 

°4,2 

(0) 

= 0, 


(0) 

= 0, 




6 (1) 

4,1 

(Z?) 

= 0, 




6 (1) 

°4,2 

09) 

= 0, 


(P) 

< 

,/3 

^sir 

■*(« 

5 

ft (2) 

4,1 

(0) 

= - 

3 j> 

Jt: 4,l 

4 

sin ( ^ ) sin 

(P), 

b (2) 

°4,2 

(0) 

= 4^ sin 2 (f) , 

b[ 2 l (0) 

= 0, 




fe (2) 

4,1 

(Z?) 

= 0, 




fe (2) 

°4,2 

09) 

= iy sin 2 (f) , 

(P) 


sin {P) 

V54,0 (0) = 

= 1, 

<^4,1 

(0) = 

= 0, 

^4,2 

(0) 

= 1, 


¥>4,3 (0) 

= 0, 


^4,4 

(0) = 

1, 

¥>4,0 (/3) = 

= 1, 

¥>4,1 

(/5) = 

= sin (/3), 

^4,2 

C/3) 

= cos (/3), 

¥>4,3 (P) 

= sin (2/3), 

^4,4 

09) = 

cos (2/3), 

V’ilo (°) = 

= 0, 


(0) = 

= 1, 

Fi.l 

(0) 

= 0, 


¥>^3 (°) 

= 2, 


¥*i|l 

(0) = 

0, 

V’im (Z 3 ) = 

= 0, 

VaI 

(P) = 

= cos (/3), 


(P) 

= - 

sin (/3), 

¥>4^3 (P) 

= 2 cos (2/3), 

¥*i!l 

09) = 

-2 sin (2/3), 

^ (0) = 

= 0, 


(0) = 

= 0, 

¥>i 2 2 

(0) 

= - 

1, 

¥>i 2 l (0) 

= 0, 



(0) = 

-4, 

¥>4 2 o (/3) = 

= 0, 

¥>i 2 l 

(0) = 

= - sin (/9) , 

¥>i 2 2 

08) 

= - 

cos (/3), 

¥>i 2 3 (P) 

= -4 sin (2/3) 

> ¥*4 2 1 

09) = 

—4 cos (2/3). 
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Substituting for n = 4 and a = 0 the derivatives above into identities (10) and (11), one obtains the transformation 
matrix 


t* ■ 

°i,3 


4, 4 

2=0, j =0 


0 \ tan (f) sin (^)~ § cos (/3) tan (§) sin (/?) 

1 1 2+4to°2(^) cos ^ + | sin (P) tan (§) cos(/3) 


0 tan ( ^ 


6sm( //g^ sin (2/3) - cos (2/3) tan (§) sin (2,3) 

4 cos^I ^ ) \ / 

cos ( 2,fl ) + sin (2/3) tan (^) cos (2/3) 


(34) 


based on which Fig. 3 shows control net configurations for the exact description of patches of some integral and 
rational trigonometric surfaces. 



Fig. 3: Control point configurations for the exact description of patches of (a) a special integral variant of 
Alfred Gray’s non-orientable Klein Bottle and of (6) a rational ring Dupin cyclide. Besides the trigonometric 
basis transformation (34), in cases (a) and (6) Theorem 2.3 and Remark 2.4 (i.e., the extension of Algorithm 
2.1) were applied, respectively. 


3.2 Hyperbolic polynomials 


Now, let a = 0 and p > 0 be fixed parameters. Using hyperbolic sine and cosine functions in expressions 
(28)-(30) instead of the trigonometric ones, we obtain the vector space of hyperbolic polynomials of order at 
most n (or degree 2 n) the unique normalized B-basis of which was introduced in [32]. In this case 

<P2n,0 ( U ) = °> i = 1, 2, . . . , 7T,, 

b-2n,0 (0) = b 2 n,2n {P) = 1, b 2n ,i (0) = b 2n ,2n-i (P) — 0, i = 1, 2, . . . , n 

and 

jn m i = / °’ j ( mod2 ) = °> <P) rn\ = ! ij . sinh (*^) > J ( mod2 ) = o, 

2n,2i-i 1 { i J ", j (mod2) =* 1, 2n,2»-i v { f cosh (ip), j (mod2) = 1, 

Jj ) mi = / ( mod2 ) = , n tt) ta\ = I iJ cos h {ip ), j (mod 2) = 0, 

F 2 n,2i K'J) \ 0j j ( mod 2) = 1, ’ ^ 2 n,2i [p) \ i> Sinh (ip ), j (mod 2) = 1 

for all i,j € {1,2,..., n}, while the higher order derivatives (f4E ( d ) > &2n i (^)}?=o ^‘=1 are specified by the 
next theorem. 


Theorem 3.2 (Hyperbolic endpoint derivatives) For arbitrary derivative order j = 1,2,... ,n, one has that 

n—r—1 r 


^2n,2r+l ( d ) 
c 2n,2r+l 


Q®< 


bt E E(-d 

k =0 1=0 


_ 1 ^n+((n — r—1) mod 2) + (r mod 2) — k—i— 1 ^2(n —i—1) + 1^ ^2r+l^ __ 


k 


^(n — k — 2r + t — 1 y — (n — k — cosh ^(2 (n — k — r) — 1) , j (mod 2) = 0, 


n—r—l r 


1 y ^ y ^ ^_-^n+((n — r— 1) mod 2) + (r mod 2) — k—t ^2(n — r— l) + 1^2r+l^ 

k =0 ^=0 

• ^(n — k — 2r + £ — 1 ) J — (n — k — sinh ^(2 (n — k — r) — 1) ^ , 


j (mod 2) = 1, 
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for all r = 0,1,..., n — 1, while 


^2n,2r ( 0 ) 

°2n,2r 

/2(n—r)\/2r\ r+(™od2)-l 

1 J2 (-l)H-(rmoda)-i ( 2 r )2J - (r 

2 r \ (n —r) + ((n —r) mod 2) —1 

V (_l)("-*0+((™-r)m°d2)-fc (2(«- r) ) ( n - r - k) j cosh {(n-r - k) 0) 

k =0 

(n—r) + ((n — r) mod 2) — 1 r+(rmod2) — 1 


= 


(V) 

" 2 2 "- 1 


1 E E ^_ -j^n+((n —r) mod 2 ) + (r mod 2) — k — £ ^2(n — r)^ ^2r^ _ 


- 2 2n-l 


£=0 


^(n — k — 3r +- 2£) J +• (n +- r — k — 2£)^ cosh ((n — r — k) (3 ), 


j (mod 2) = 0, 


(2{n-r)\ (2r\ (2(n-r)\ r +( r mod 2 )~ 1 

{ n -W Kr) + l 2 L~-V E (-l) r+(rmod2) -'( 2 ;)2^(r-^sinh(2(r-^)/3) 

£=0 

(n — r) + ((n — r) mod 2) — 1 r+(rmod2) — 1 

E E ^_-j^n+((n — r) mod 2) + (r mod 2) — k — £ ^2(n —r)^ ^2r^ 

fc=0 £=0 


2 2n-l 

• ^(n +- r — k — 2£) J — (n — k — Sr +- 2Z) J J sinh (2/3 (r — ^)) , 
for all r = 1,2,..., n and i (/3) = (-1) J 6 ^ 2n _ i (0) , z = 0,1, ..., 2n, j = 0,1, ..., n. 


j (mod 2) = 1 


Example 3.2 (Second order hyperbolic polynomials) Using hyperbolic sine, cosine and tangent functions in¬ 
stead of the trigonometric ones that appear in Example 3.1 and applying the second order hyperbolic normalized 
B-basis [32] with the shape parameter /3 > 0, one can easily construct the hyperbolic counterpart 



4, 4 

t=0, j =0 


1 1 

0 itanh(f) 

1 1 

0 tanh (4) 

1 1 


1 


1 


1 


3 sinh(/3) 
2+4 cosh 2 ^ ^-) 

3(l+cosh(/3)) 
2+4 cosh 2 ^ ^) 

6 sinh(/3) 
2+4 cosh 2 ( ^-) 

6 cosh(/3) 
2+4 cosh 2 ^ ^) 


sinh (j8) — | cosh (/3) tanh ( [ ^ sinh (/3) 
cosh (/3) — I sinh (/3) tanh ( + cosh (/3) 
sinh (2/3) — cosh (2/3) tanh (^ sinh (2/3) 
cosh (2/3) — sinh (2/3) tanh (^ cosh (2/3) 


(35) 


o/ the trigonometric basis transformation ( 34 ), the structurally difference of which consists in the highlighted operators. 


3.3 A class of mixed spaces 


In order to be as self-contained as possible, we recall the construction process [6] of the normalized B-bases for 
a family of mixed EC vector spaces of functions. 

Let a = 0 and /? > 0 be fixed parameters and consider the homogeneous linear differential equation 

n+1 

y: 7 i v ^ ( U ) = 0, 7i G M, ue [0,/3] (36) 

2 = 0 


of order n -(- 1 with constant coefficients and assume that its characteristic polynomial p n +i (r), r £ C is an 
either even or odd function such that r = 0 is one of its (presumably higher order) zeros. Hereafter we assume 
that the ordinary basis (1) corresponds to the system of those linearly independent functions that are implied 
by all (higher order) zeros of p n + i, i.e., is the (n + l)-dimensional vector space of functions that is formed 
by all solutions of (36). Under these conditions, 1 € §S+, moreover the space is also invariant under 
reflections and consequently under translations as well, i.e., for any function / £ §5+ and fixed scalar t £ ]R 
the functions g T ( u ) := / (r — u) and h T (it) := / (u — r) also belong to §° ,/3 . 

Following [6], one can both to determine (or at least to numerically approximate) the range of the shape 
parameter (3 > 0 for which is an EC space and to construct its normalized B-basis as follows. Denote by 


Wr 





n, n 

2 = 0 , j =0 


the Wronskian matrix of those particular integrals 


n 

v n,i ’ = y ^ Pi,j^Pn,j ^ ^>n ^ i ^ — 0, 1, . . . , 71 
3 =0 


(37) 


(38) 
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of (36) that correspond to the initial conditions 

f v® (0) = 0, j = 

\ (0) = 1, (39) 

{ v n.'} (P) = °. j = 0,...,n-l-i, 


i.e., the system {v Ut i («) : u £ [0, /3]}”_ 0 is a bicanonical basis on the interval [0,/?] such that the Wronksian (37) 
at u = 0 is a lower triangular matrix with positive (unit) diagonal entries. 

Consider the functions (or Wronskian determinants) 


w n,i (u) := detw [ Vnii>Vnti+1 ,..., Vn<rl ] 0) > 

define the critical length 



/3n : = min min{|«| : w n> j (u) = 0, u ^ 0} 

i =LfJ-LfJ+ 1 ’-’" 


(40) 


(41) 


and, in what follows, assume that /3 £ (0,Pn) is an arbitrarily fixed shape parameter (we write pn = +oo 
whenever the Wronskian determinants (40) do not have non-zero real zeros). Under these conditions, is a 
reflection and translation invariant EC space that also has a unique normalized B-basis, since 1 £ 

Consider the Wronskian matrix W^ v v v q j (/ 3 ) of the reverse ordered system {v n ^ n -i (u) : u, £ 

[0, /3]}™_ 0 at the parameter value u = /3 and obtain its Doolittle factorization 


L-U = Wr 




03 ): 


where L is a lower triangular matrix with unit diagonal, while U is a non-singular upper triangular matrix. 
Calculate the inverse matrices 


U- 1 :^ 


Po.o Pop ■ 

’ MO ,n 


Ao.o 

0 • 

■ 0 

0 pip ■ 

’ Ml ,n 

, L~ x := 

Al,0 

Ai,i • 

■ 0 

0 0 ■ 

■ Mn,n. 


_ An,0 

An,l ‘ 

■ An,n_ 


and construct the reflection invariant normalized B-basis 


K>°,P _ 
Un. — 


defined by 
and 


^ bn,i (r^) — i,O b n,i ('a) ■ U £ [0, /3] ^ 

[bn,n («) b n ,n -1 («)••• b 0 (it)] := [v n ,n («) l>n,n-l (u) • ■ • V n ,0 («) ] ' U~ 


(42) 


[A 0 ,o Ai,o A„, 0 ] ~L~ X - [10---0] . 

Since the EC space is invariant under reflections, one has that 

bn , i (n) — b n n —i (/I It) , i — 0, 1, ... , — J : 

i.e., we only need to determine the half of the basis functions (42). 

{ / -\ / n», n 

<Pn,c (0) , <Pn,c (P) f are already 

’ ’ J c= 0, j =0 

known, one has to substitute the parameter values u = 0 and u = p into the derivative formulas 

= . (y\ 

n,n—i V / b u n,n—i V / 
i 

= K,0 h'r^n-r («) 


i n 

— ^i,0 ^ ^ Mr ,i ^ ] Pn—r,cPn^c ( u ) ? ^ — 0, 1, , |^ —J , 


(43) 


r =0 c=0 


b n] ( w ) = (- 1 )- 7 b ni-i (P-u), i = o, 1,..., UJ (44) 

in order to determine the entries (10) and (11) of the transformation matrix that maps the normalized B-basis ()2) 
of the (generally mixed) EC space E>n^ to its ordinary basis. 

From computational and algorithmic viewpoints, formulas (43)-(44) are significantly easier to both evaluate 
and implement for fixed values of the shape parameter /3 than to calculate the general determinant based 
formulas (12)-(14) or to differentiate the integral representation described e.g. in the special case [20] and 
references therein. 
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Example 3.3 (EC spaces generated by characteristic polynomials) Consider the characteristic polynomials 


Pn+1 ( r ) = r n+1 , p (n+ 1)2 (r) = r ” +1 J} (r 2 + w t) , P {n+1 f (r) = r n+1 J} (V 2 - 


n+1 —fc 


fc=l 


fc=l 


where parameters {<^fc})j =1 are pairwise distinct non-zero real numbers. In these cases, one has that for appropriately 
selected definition domains [0, (3\ the vector spaces 

P° ,/3 := ({l,u, • • • ,u n : u G [0, /3]}) , dimP°’ /J = n + 1, 

AT n(»+ 2 ) :=Frl u cos(u k u) ,t/sin(w fe u) : u G [0,£] , dimAT°’ ( ^ +2) = (n + l) 2 , 

and 


AE n(n+ 2 ) :=P„U ^ ji/ cosh (ajfcit) , u sinh (u k u) : u G [0, , dim AH°' ( ^ +2 ) = (n + l) 2 , 

respectively, are reflection invariant EC spaces that also possess unique normalized B-basis functions. As special cases, 
the vector spaces of trigonometric and hyperbolic polynomials of order at most n correspond to the characteristic 
polynomials 

n n 

P 2 n +1 (r) = r (r 2 + k 2 ^j and P 2 n+i (r) = r JJ (r 2 - k 2 ^j , 
k =1 k =1 

respectively, i.e., uif. = k for all k = 1,2,... ,n. However, in these two latter cases it is much easier to apply Theorems 
3.1 and 3.2, respectively, than to evaluate the required endpoint derivatives by means of formulas (43)-(44)- Concerning 
the characteristic polynomial p n +l (r) = r n+ , Example 3-4 and Appendix A provide further details. 

Example 3.4 (Traditional polynomials) The system < 6 r(? ; (u) = (+)vf (1 — u) n ~ l :«£ [0,1]| of Bernstein 

polynomials of degree n is the normalized B-basis of the EC space (u) = u 1 : u G [0,1]. In this case one 

has that 



*!, j = i, 

0, j ^ *, 

0, 



(-i y-'-ji ■$)■({), 


a 

0, 


, i > j > 0, 

i<j<n, 


i > j > 0, 
i< j <n. 


&n?i (0 


(-lMLi(°)> i = 0,1,... ,n 


(45) 

(46) 


for all j = 0,1,... ,n. As it is proved in Appendix A, the substitution for a = 0 and /3 = 1 of these derivatives into 
formulas (10) and (11) leads to the expected closed form of the classical transformation matrix of entries 




( j ) 

7sy, j = + n, 

0, j = 0,1,... ,i — 1, 


(47) 


where i = 0,1,..., n. 

Naturally, characteristic polynomials may also have (conjugate) complex roots of higher order multiplicity 
and with non-vanishing real and imaginary parts, which may lead to mixed (algebraic) exponential trigonometric 
EC spaces as it is illustrated in the next example. 

Example 3.5 (A 5-dimensional exponential trigonometric space) Let uj > 0 be a fixed parameter and con¬ 
sider the 5 th order homogeneous linear differential equation 

v ^ (u) — 2 ^u; 2 — 1^ v^ (u) + (uj 2 + 1 j t/ 1 ) ( u ) = 0 (48) 

with the odd characteristic polynomial 

P 5 (r) = r 5 — 2 (uj 2 — 1^ r 3 + (uj 2 + 1 j r 

= r(r - ( uj - i)) (r - (-uj + i)) (r - ( uj — i)) (r - (w + i)), 
where i = y—l. It follows that the vector space formed by all solutions of (48) can be spanned by the ordinary basis 

-pf = jipo (w) = 1, <pi (u) = e~ ulu cos (u), g >2 (u) = e~ ulu sin (u), (49) 

ifi 3 (u) = e““ cos (u), <fi 4 ( u ) = e uu sin (u) : u G [0,/?]} , /3 G (O,/^) , 

where denotes the corresponding special case of the critical length (41). In order to avoid lengthy cumbersome 
formulations, in this case we provide only a numerical example the values of which can be verified by means of 
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Listings B.l and B.2 of Appendix B. Assume that the growth rate lo = ^ and the shape parameter /3 = < flX are 

fixed. If one intends e.g. to represent the arc 


e ulu cos (u) 
esin (u) 


1 

0 


</>3 («) + 


1 ¥>4 (u) , « G [o, 0] 


(50) 


of a logarithmic spiral by means of the normalized B-basis { 64,4 (u) : u G [O,0]}^_ o 0 / the underlying reflection in¬ 
variant EC space S 4 " 8 , then one has to construct the system (42) as follows: 

— at first, one has to determine the transformation matrix 


r 14,4 

[Pi,j\i=0, j—0 


- 0.8038 
-1.2028 
1.4484 
-1.1191 
0.9779 


-0.7765 

4.2007 

-4.8494 

2.8895 

-0.4889 


-2.4061 

2.2514 

0.5805 

-2.3598 

2.2781 


0.9728 

-2.9979 

3.4009 

-1.7704 

-0.4889 


1.0761" 

-0.4876 

-1.4558 

2.8542 

-2.2781. 


that maps the ordinary basis {^ 4 ,?: (u) : u G [0, /3 ]}^_ 0 to the particular integrals { 14 ^ (u) : u G [0,/3]}f =0 of the 
form (38) which fulfill the initial conditions (39); 

— next., one has to obtain the Doolittle LU-decomposition of the Wronskian matrix 


W, 


[^4, 


4.>'4,3>V4,2,'W4,l>t'4,o 


]( 0 ) = 


L-U = 


1.2166 -0.0000 -0.0000 0 -0.0000 
1.3923 -0.9304 -0.0000 0 -0.0000 

0.6629 -1.6472 1.0000 -0.0000 0 

-0.8886 -0.5383 1.7705 -1.0748 0.0000 

-1.5551 2.1130 -0.4963 -1.2300 0.8219 


of the reversed ordered system {f 4 , 4 -i ( u ) :«G [0, at u = /3, i.e., 


■ 1.0000 

0 

0 

0 

0" 


"1.2166 

-0.0000 

-0.0000 

0 

-0.0000" 

1.1444 

1.0000 

0 

0 

0 


0 

-0.9304 

0.0000 

0 

-0.0000 

0.5449 

1.7705 

1.0000 

0 

0 

,U = 

0 

0 

1.0000 

-0.0000 

0.0000 

-0.7303 

0.5786 

1.7705 1.0000 

0 


0 

0 

0 

-1.0748 

-0.0000 

_-1.2782 

-2.2711 

-0.4963 1.1444 1.0000. 


0 

0 

0 

0 

0.8219. 


— then, by using the essential parts of the inverse matrices 



■ 1.0000 = Ao.o 0 0 0 0" 


’0.8219 = /x 0,0 —0.0000 = /xo,i 0.0000 = /xo ,2 * *" 


-1,1444 = Ai, 0 * 0 0 0 


0 —1.0748 = /xi ? i 0.0000 = /xi t 2 * * 

L- 1 = 

1.4812 = A2 ,o * * 0 0 
★ ★ ★ ★ 0 

, u~ l = 

0 0 1.0000 = /x 2,2 * * 

0 0 0 ★ ★ 

0 0 0 0 *. 


and the reflection invariant property of the vector space, one has that 


64,4 (u) = Ao,oMo,ot> 4,4 (u), 

&4,3 («) = Ai,o (/to,1^4,4 (u) + /(l,lt>4,3 («)) , 

64,2 (u) = A 2,0 (/t0,2f4,4 (u) + /(l,2f4,3 («) + /t2,2f4,2 («)) , 


64,0 (u) = 64,4 (fi —u), 
b 4 ,i (u) = 64,3 (f3 - u) , (51) 

u G [0,0] . 


Fig. 4(a) shows the image of these normalized B-basis functions. 



Fig. 4: (a) Exponential trigonometric normalized B-basis functions of order 4 which correspond to the shape 
parameter (3 = ^ and growth rate ui = ( b ) Control point based exact description of different arcs of the 

logarithmic spiral (50). 
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Using the higher order derivatives of the ordinary basis functions (49) at u = 0 and u = 13, one can also easily 
evaluate the higher order derivatives of the obtained normalized B-basis functions by means of formulas (43)-(44)- 
Substituting these derivatives into (10)-(11), one also obtains the transformation matrix 




4, 4 

i= 0, j=0 


■ 1.0000 1.0000 1.0000 1.0000 1 . 0000 " 
1.0000 0.9073 0.2057 -0.3859 -0.6560 
0 0.8738 1.0520 0.9871 0.3787 

1.0000 1.0927 0.4593 -0.4605 -1.1433 
0 0.8738 1.3386 1.5980 0.6601. 


that is required for the control point based exact description (17) of any arc of the logarithmic spiral (50) that is defined 
over an interval of length /3 (see Fig. 4(b))- Observe that from algorithmic and implementation viewpoints, the steps 
above are much easier and more efficient to perform than the evaluation of other possible integral or determinant based 
representations. As long as parameters w or ft are not modified, the calculations above do not have to be reevaluated. 

Example 3.6 (Quadratic algebraic trigonometric functions) The normalized B-basis 

= { 64,0 (u) = b 4>4 (ft — u ), 64,1 (u) = 64,3 (ft — u), 

b 4 ,2 (m) = cf 2 ^2/3 (sin (u) — sin (ft)) — 2/? (1 — cos (ft)) u + ft 2 + 2/3sin (/3 — u) — /3 2 cos (ft — u) + 

+/? 2 (cos (ft) — cos (u)) + 2 (1 — cos (ft)) u 2 + ft (ft — u) u sin (ft)) , 

b 4 ,3 ( u ) = cf 3 (2 (ft — u) + 2 (sin (u) — sin (ft)) + 2 (u cos (ft) — ft cos (u)) + 2 sin (ft — u) + 

+ft 2 (u - sin (u)) - (ft - sin (ft)) u 2 ^ , 

64,4 (u) = cf j4 (2 cos (u) + u 2 - 2^ : u G [0,/3] j , ft G (0,ft 4 ) , ft 4 = 2n 

of the EC space 8 °’^ = (J r 4 ’' 9 } = ({ 1 / 34,0 (u) = 1, 1 / 34,1 (u) = u, 1 / 34,2 (u) = u 2 , 1 / 34,3 (u) = sin (u ), 1 / 34,4 (u) = cos (u) : 

u G [0, /?]}) of algebraic trigonometric functions can also be constructed, e.g. by using either the differential equation 
based iterative integral representation published in [20] and references therein or the determinant based formulas of 
[22, Theorem 3-4]- The critical length ft 4 = 27r was determined in [6, Section 5] or [8, Proposition 3], while positive 
scalars 

fj 4 — 4 cos {ft) — 2ft sin {ft) 

4.2 {ft 2 — 4 cos {ft) — 4ft sin (ft) + ft 2 cos (ft) + 4) 2 ’ 

ft = _ 2 (ft - sin (ft)) _ 

4.3 (2 cos (ft) + ft 2 — 2) (ft 2 — 4 cos (ft) — 4ft sin (ft) + ft 2 cos (ft) + 4) ’ 

c P = 1 

4 ’ 4 2 cos (ft) + ft 2 - 2 

are normalizing coefficients. Applying Theorem 2.1 with the settings above, one obtains the transformation matrix 



4, 4 

i= 0, j= 0 


1 1 


1 


1 


1 


r. 2cos/3+/3 2 — 2 
2(/3 — sin /3) 

0 0 

n 2cos/3+/3 2 —2 
U 2(/i—sin ft) 

1 1 


(2—2 cos [3 — [3 sin [3)[3 
4—4 cos [3 — 2/3 sin (3 
f3 2 —4 cos [3—4/3 sin [3-\-(3 2 cos (3+4 
2 — 2 cos [3 — /3 sin [3 
(2—2 cos /3 — )3 sin /3)/3 
4—4 cos /3 — 2/3 sin [3 
(2 sin /3 — [3 — [3 cos /3)/3 
4—4 cos [3 —2/3 sin /3 


-o _ 2 cos /3-\-/3 2 —2 

P 2(/3—sin/3) 

„2 _ (2 cos /3+/3 2 —2)/3 
^ (/3-sin/3) 

(2 cos /3+/3 2 —2) sin(/3) 
COS(,p/+ 2(/3—sin f)) 


ft 
ft 2 

sin (/3) 
cos (/3) 


(52) 


that maps to J-, based on which Fig. 5 illustrates the control point based exact description of cycloids and 
helices of different shape parameters, while Fig. 6(a) shows the control net of a cylindrical helicoid. 


Remark 3.1 (Hybrid EC B-surfaces) Naturally, one can also combine different types of normalized B-basis func¬ 
tions in order to describe hybrid surfaces as it is shown in Fig. 6(b) that illustrates the control point based exact 
description of a hyperboloidal patch. 


4 Proof of main results 

Proof of Theorem 2.1. The linear transformation J=0 that maps the normalized B-basis 6“’^ of the 

vector space §“’ /3 to its ordinary basis will be constructed by mathematical induction on the column index 
j or n — j, where j = 0,1,..., (5J ■ Using one of the properties (5)—(8) , at each step j we will compare the left 
and right side of the jth order derivative of the matrix equality (9), thus obtaining an iterative process that is 
outlined in Fig. 7. 

First of all, observe that 

to,j = 1, Vj = 0,1,... ,n 

and 

li,obn,i (^) t/,o Tn,i (^) 3 li.nbn.n (ft) ti,n Tn,i (ft) 3 l 0, 1, . . . , Tl, 
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Fig. 5: (a) Cycloids and ( b ) helices of different shape parameters described by means of Theorem 2.2 and of 
the basis transformation detailed in Example 3.6. 




(■ b) 


Fig. 6: Control point based exact description of (a) helicoidal and of ( b ) hyperboloidal patches, respectively. 
The control net of the patch (a) was constructed in both direction by means of basis transformations of the 
type (52) with different shape parameters. In case of patch (6) the control point configuration was obtained by 
using both the trigonometric and hyberbolic basis transformations (34) and (35), respectively. 



property (5), i = 0, j = 0, 1, . . - , n, « = 8 
| property (6), j = 0, i = 0, 1, . . . , n 
3 properties (7)-(8) for j — 1, i = 1, 2, . . . , i 

| properties (7)—(8) for j — 2, i — 1, 2, , . . , i 

| properties (7)—(8) for j = 3, i = 1, 2, . . . , i 

| properties (7)-(8) forj = 4= |_^J> * = 1 > 


Fig. 7: Outline of the proof. 
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due to the partition of unity property (5) and to the endpoint interpolation property ( 6 ), respectively. Using 
forward substitutions, the elements of the columns [i" 7 ]" =1 , j = 1 , 2 ,..., [§J are iteratively determined by 
differentiating the matrix equality (9) with gradually increasing order and applying the Hermite conditions (7) 
at u — a. In order to formulate a mathematical induction hypothesis, let us consider some special cases. When 
j = 1 one obtains that 

fin} («) = 4>5 (a) + ( a ) > * = °, 1 , ■ • ■, n, 

where b^ \ (a) 0 and for the special subcase * = 0 one has that 

bn,o ( Q ) + b n\ («) = fin}o («) = 0, 


b nl («) = - 6^1 (a) ■ 


4-' L _ _ 

M “ 1,(1) 


1 


&5 («) 

i = 1,2,..., n. 
For j = 2, we have that 


“ f £o 6 5 («)) = m 1 (</>5 («) + <p n ,i (a) (a)) 

V ' b n,l (“) 


( 1 i 

- (a) + (1) . 

6 „,1 («) 


'fin} (a) = tifib { n}o (a) + tllbnl (a) + t” 2 &5 («) 
where b^\ (a) yfi 0 and for the special subcase i = 0 we obtain that 

b n}o («) + b n}i (a) + bfy (a) = <p { *} 0 (a) = 0 , 
b n}o («) + b n} («) = -&i 2 2 (a), 


i.e., 


+ '. L _ 


- 1 /U( 2 ) 


6 5 («) 


(^5 («) - («) - f”i^ 2) i («)) 


iy— («) - <p n ,i (a) &5 (a) - fip n ,i (a) + ^ \ ^ b n\ (a)] 

, 2 ( a ) \ \ Vi(a)/ J 


fc(2) : W V 

(1) , x (2) 

= fin,i (a) - ^ ^ &5 (a) + 5 ’C ; ; , i = 1 , 2 ,..., n. 


& i 2 2(a) ^i(«) 71,1 &i 2 2(a) 

In case of j = 3 one obtains that 

fin} («) = tifibn ,o (a) + (a) + *" 2 65 («) + b h b n}z (a), * = 0 , 1 ,..., n, 

where &5 (a) 7 ! 0 and for the special subcase i = 0 , one has that 

b n}o («) + 6 5 («) + b n} 2 («) + 6 5 («) = finfi («) = 0, 


~ b n}i («) = b n}o («) + 6 i 3) i («) + b n}i («). 


i n — (,S 3) 

l ’ 3 ^(a) 


{fin} (a) - tiA% (a) - t"i&5i ( a ) - ^2^5 (“)) 

“ - 1 - ^ 5 3 ] ( Q ) " W* («) b n} 0 (a) - (Vn,i (a) + ^ ! x ) 6.5 ( Q ) 

V b n, 1 (a) / 

*>£!(«) .(2) , X . ^(“A t (3) f A 

y (a) 

6 $ (a) V 6 5(a) V b n ,2 (“) 


fe i 3 3 («) 

- Vn,i (a) - 


n,t v-v | , (3) / x b 5(a) ^5 (a) ^ l 5w,(3) , fin} ( Q ) 




(a) + 


, ( 2 ) , S n, 2 \ J ' ( 3 ) ■ 

6)1,2 (a) / 6)J (a) 
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One can observe that expressions corresponding to these special cases are in accordance with formula (10). 
Now, fix the column index j £ {l, 2 ,..., |_§J — 1} and assume that formula ( 10 ) is valid up to the selected index 
j and we will also prove it for j + 1. We can proceed as follows: 


3 + 1 

Vnf 1 ' 1 ( a ) = E *£ 7 & «,7 ^ («) > « = 0, 1, . . . . 

7=0 


where 6 ^ 1 "^ (a) 0 and for the special subcase i = 0 we have that 


i+i 


E ( Q ) = ‘fin,o ' 1) ( Q ) = 

7=0 


i.e., 


~ b nJ + l («) = E 6 «.7 1) («) > 
7=0 


k ’ j+1 ~ 6 (j+i) 

n,j + l 


(a) 


<r («)-E*M£r w 


7=0 


„o+i) 


(«) 


b ( ' j+1) (a) b (j+1) (a) ^ 

U n,j +1 V a / u n ,7 + 1 V a / 7 =l 


E &n ,7 (“)*£• 


n,j+l v^v 7=0 


^ +1 > 
^ n,2 


(a) 


b (i+1) (a) (a) ^ 

U n,j +1 VW u n.7 + l 7= 


E 1} (a) I <Pn,i (a) + 


n,j +1 


7=0 


,(t) 

In) 


fin i («) 1 


7 - 1 ,.W 
n.i 


bn ,7 (a) (<*) J^l (a) 




7—7—1 

+ E (-1)' E 

£=1 r</c 7 ) i<fc 7 ) 2 <...<fe 7) ^<7 


, (fc 7 ,l) / \ 2 ) / \ » ( fc T.3) 


(fc 7 ,r) 


«?*' («> ‘is;; <«> ck <«>•••c;;;.. <«> ‘S,, <«) 


j (^7,3.) / \ , (^ 7 , 2 ) / \ / \ 

0 , (a)0 , (a)... b 1 (a) 

n,fc 7) 1 V / n,/c 7) 2 ' ' n,k-y y £ V ' 


W , 1 

{,W+!) ( Q ) Q ftC + 1 ) (q) 

°ra,j+l W °nj+l 7=0 

7 — 7—1 

+ £ (-o' 


6^7 («) 6^7 (a) ^ b$ r (a) V 


<=1 


E 

7’</C7 ; 1 < kry y 2 < • • • < ^7 ,£ <7 


0 ’ 7 ^ V ^7 ( a ) ^7 ( a ) ^ 

& a) («) («) («)•■■ («) («)' 
\ (u \ (u .\ 


' ° 1 ,v 7,1 t- ; < V7 , ^ _ 1 t, t ,1 '7,1 — 1 

7 (^7,i) / \ 7 (^7,2) / \ 1,(^7.^) / 

o v , (a) o , (m ... o v , (o 

n,/c 7) i v / n,/c 7) 2 ' ' n,/c 7) ^ v 


1 


*£! ( q ) /" . 0 + 1 ) 


^ , ' 0+1) («) S (a) 


J-r 


"J+l 


^ (a) 6 ^ (a) &J& (a)... ^ (“) («) 


(«) + 




hO+i) 


+E(- 1 ) E 

i=l r<fci<fc 2 <..<fe£<) + l 


‘Si 1 .(«)<!,(“) • 15W 


.(b) 


(fcr) 


+ 


r n,l 


(a) 


E;E'0 


which means that our induction hypothesis is correct for all column indices j = 1 , 2 ,..., [_5J- 

Using a similar technique based on backward substitutions, the entries of the columns [t" n _j]” =1 , j = 
1,2,..., L^J can also iteratively be determined by differentiating the matrix equality (9) with gradually in¬ 
creasing order and applying the Hermite conditions ( 8 ) at u — j3. After correct reformulations one obtains 
exactly the formula ( 11 ). ■ 


Proof of Theorem 2.2. Using Theorem 2.1 and Corollary 2.1, the fth coordinate function (£ = 1, 2,..., 8) of 
the ordinary integral curve (16) can be rewritten into 

n n 

C ( u ) = E^n,i ( U) = E P3 hn d ( U ) ’ Vm £ [ a ’^] > 

2=0 j =0 

n 

Pj = E i = 0 , i,..., rr. 

2=0 


where 
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Repeating this transformation along all coordinate functions and collecting the coefficients of the normalized 
B-basis functions, one obtains the vertices p y = of the required control polygon. ■ 

Proof of Theorem 2.3. By means of Theorem 2.1, one can construct for all r = 1,2 the regular transformation 

matrix \t nr that maps the normalized B-basis Bn r ’^ r of the vector space to its ordinary 

L lr ’ 0r ii r =o, j r =o 

basis J r nl’^ r ■ Observe that the £th coordinate function (£ = 1, 2, 3) of the ordinary integral surface (18) can be 
written in the form 


(T£ 2 


<Jt 2 


* ( U ) = E II I X ifvn r ,ir- K) I = E II E P f jl Cb ^r,:ir K) 

C=lr=l \i r =0 / C=l r =l \jr=0 ) 

ni n 2 / o-f 2 \ n i n 2 

= El El I El Jl[ Pjf I ^ni,ji ( u l)bn 2 ,j 2 ( w z) = El El Pji,j 2 ^ ) n 1 ,j 1 ( u l) b n2 j 2 ( 112 ) 

.71=0.72=0 yc=lr=l J j 1= 0j 2 =0 

for all u = [u r ] 2 r=1 G [«i,/3i] x [ 02 ,^ 2 ], where 


<7 £ 2 


=e n 4 C 


and the values 


C=1 r=l 


= E EEb,-’ C = !, 2,..., at, jr = 0,1,..., n r , r = 1,2 


can be obtained by means of Corollary 2.1. Repeating this reformulation for all coordinate functions and 
collecting the coefficients of the product of normalized B-basis functions, one obtains all coordinates of all 

control net points Pj lt j 2 = |pj ljj2 J . ■ 

Proof of Theorem 2.4. Let i £ {1,2,..., n}, j € {1,2, J}, r G {1, 2,... ,j — 1} and l G {1, 2,... ,j — r — 1} 
be fixed indices at the moment and consider formula (10). There are O - ), -1 ) pairwise distinct strictly increasing 
sequences r < k\ < k 2 < ■ ■ ■ < kg < j of length l between r and j that can also be stored in permanent lookup 
tables, since they are independent of the applied normalized B-basis. In case of each of these sequences one has 
to evaluate the fraction 


(a) b^ 2 k [ (a) b™ (a) ■ ■ ■ ■ ■ ^-1 (a) (a) 

6 (fcl) fad 6 (fe2) (a) • ■ b (ke) (a) 

U n.k , W "ni, W ' ' ' "n.fa W 


*,( fe 2) 


u( fc 3) 


Of) 


uti) 




that includes 21 flops (i.e., I multiplications in the nominator, £—1 multiplications in the denominator and 1 
division). Thus, the total number of flops required for the evaluation of the summation 


s j,r,i E 

r<k 1 <k 2 <...<ke<j 


(53) 


equals 


2£• 


j - r - 1 


+ 


j-r-1 


- 1 


for each fixed values of £ = 1,2, ...,r — 1, where the last term in the parentheses appears due to additions 
that have to be performed in (53). If one considers all possible values of £ and observes that (— l) e is just an 
alternating sign (implying either addition or subtraction), the number of flops performed during the evaluation 
of the expression 

j—r—l 

bnjr («) + E (-l/ Sj, r ,£ 


9j,r — 


bn,r («) 


7—7 - 1 

E (*■ 

f=\ 


j - r - 1 


+ 


j-r-1 


- 1 + [1 + O' - r ~ 2)] + 1 = 2 J ' 


0 - r ) 


that consists of the evaluation cost of all { s j,r.e}j = i , of 1 + (j — r — 2) = j — r — 1 additions and of 1 division. 

Observe that values { g^ r are independent of the row index i for all fixed column indices j = 1,2,..., [§J“ 
i.e., they can be evaluated and stored in a temporary lookup table by performing 


3 -1 


E 2 i_r_1 '0 - r) = 2 j ~ 1 ■ j - 2 j + 1 
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flops and later they can be reused for the evaluation of all quantities 

i -1 

„(r) 


hi,3 '■= E Vnl («) • 9j,r, i = 1, 2, ■ • ■ ,n. 


Thus, independently of i. the calculation of hj j takes an additional 2 j — 3 flops (i.e., j — 1 multiplication and 
j — 2 addition) for all fixed values of j = 1, 2,. . ., |y J. Finally, each of the column entries 


t'i.j — *Pn,i (tt) 


h i,3 - Vnl ( Q ) . 


b n) («) 


i = 1,2 ,n 


can be evaluated by means of 3 additional flops for all fixed values of j = 1, 2,..., |^J. 

Thus, the total number of flops required for the evaluation of unknown entries the general 

transformation matrix is 

E ( 2J_1 • 3 ~ 2 J ' + l) + n ■ [(2 j - 3) + 3] = 2 ^J ( |2j _ 3 ) + [£j + 3 + n ( |jj + l) . (54) 


3 =1 


3 = 1 


Since the structure of formula (11) is very similar to that of (10), one can conclude that for odd numbers n 
the total computational cost of all unknown entries of the general transformation matrix is twice of (54), while 
for even values of n the entries of the middle column do not have to be reevaluated by means of (11), i.e., in 
this latter case the partial computational cost (54) has to be increased by 


E (2 i - , .j-J< + l)+„.^[(2j-3)+3] = 2*->(!-4) + 2+2 + ^.(!-l) 


3=1 3 = 1 

which leads to the final expression (26). 


Proof of Theorem 3.1. In order to determine the higher order derivatives of normalized B-basis functions (29) 
at the endpoints of the interval [0,/l], we will make use of trigonometric identities 


sin 2r+1 (0) = 22 ?+tX!(" 1 ) 


-e 2r+l 


e=o 


• 2 r 

sin 


w = i( 2 ;) + ^E(- i r W2r 

£=0 


sin ((2 (r - i) + 1) 6 ), 
cos ((2 (r-l))9), 


where r € N and 0 € R. E.g. if % = 2r + 1 (r = 0,1,..., n — 1), then 


&2n,2r+l (m) _ gin 2(ri-?--l) + l ( P ~ u \ airi 2r+l 


v 2n,2r+l 


2 

n—r—1 


sm 


© 


dtlHl E (->)' 


22 (n— r—1 ) + 


k=0 


-r-i-fe ( 2 (n- r - 1) + l\ . 


sin I (2 (n — r — fc) — 1) 


ft — u 


22r+l 


El- 1 )' 


-£ / 2r + 1 \ . 


£=0 
n—r—1 r 


hn E E<-»> 


n+1- 


((2 (r~£) + 1)|) 


_ fc _W2 (n-r-l) + lW2r + l 


2 2(n-l) 


fc=0 £=0 


sin I (2 (n — r — fc) — 1) 


p — u\ . 


sin ((2 (r-£) + 1) 


n—r—1 r 


= 02^1 E EM) 


fc=0 £=0 


n+l-fc-f |2(n — r — l) + l^|2r + l 
fc 


cos ( (n — fc — £) u — (2 (n — r — fc) — 1) 


— cos I (n — fc — 2r + l — 1) u — (2 (n — r — fc) — 1) 


P 


from which follows that 


b U) 

u 2n,2r+l 


0) 


1 


u 2n,2r+l 


2 2 n ~i 


n —r—1 r / \ / 

y^ y ^y+i-(fc+t) /2 (n-r-l) + lj/2r+l 

a-=o r=o V ^ / V 
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(n — ( k + £)) 3 cos I (n — (k + £)) u — (2 (n — r — k) — 1) — + 

\ z z 

— (n — k — 2 r + £ — l) 3 cos ^(n — k — 2r + i — 1) u — (2 (n — r — k) — 1) ^ 

for all j > 0. Substituting u = 0 into the last expression, one obtains exactly the formula (31). If i = 2r 
(r = 0,1,..., n), then one can proceed analogously. ■ 


Proof of Theorem 3.2. In order to determine the higher order derivatives of the hyperbolic counterpart of 
the normalized B-basis functions (29) (see also Subsection 3.2) at the endpoints of the interval [0, /?], one can 
follow the steps of the proof of Theorem 3.1 by applying the hyperbolic identities 


sinh 2r+1 (9) 
sinh 2r ( 6 ) 


22r+l 


£(-D 


■+(rmod2)—/ / 2r + 1 \ . 


sinh ((2 (r — £) + 1) 9). 


. . r+(rmod2) — 1 . 

1 2r Ul V r_n r +( r mod 2 )-e I 2r 

2 2r r + 2 2r 2^ ^ > 


e=o 


cosh ((2 (r — £)) 9) 


and basic properties 


sinh ^ (uju) 

cosh ^^ (uju) 
2 cosh ( 6 * 1 ) cosh ( 6 * 2 ) 


J uj 3 sinh (uju) , j (mod 2 ) = 0 , 
\ uj 3 cosh (uju) , j (mod 2 ) = 1 , 

J uj 3 cosh (uju) , j (mod 2 ) = 0 , 
( uj 3 sinh (uju) , j (mod 2 ) = 1 , 

cosh ( 6 *i + 6 * 2 ) + cosh (9\ — 6 - 2 ) 


of the hyperbolic sine and cosine functions, where reN and 9,6 i, 02 ,w G R. ■ 


5 Final remarks 

As listed in Section 1, concerning geometric modeling, the normalized B-bases (of EC spaces that also comprise 
the constant functions) ensure many optimal shape preserving properties and algorithms. Moreover, they may 
also provide useful design or shape parameters that can arbitrarily be specified by the user or the engineer. 
In Section 3, we have seen that polynomial, trigonometric, hyperbolic or mixed EC spaces allow us to obtain 
the control point based exact description of many (rational) curves and surfaces that are important in several 
areas of applied mathematics. The investigated large classes of vector spaces also ensure the description of 
famous geometrical objects (like ellipses; epi- and hypocycloids; Lissajous curves; torus knots; foliums; rose 
curves; the witch of Agnesi; the cissoid of Diodes; Bernoulli’s lemniscate; Zhukovsky airfoil profiles; cycloids; 
hyperbolas; helices; catenaries; Archimedean and logarithmic spirals; ellipsoids; tori; hyperboloids; catenoids; 
helicoids; ring, horn and spindle Dupin cyclides; non-orientable surfaces such as Boy’s and Steiner’s surfaces 
and the Klein Bottle of Gray). 

EC bases of type (1) represent a large family of vector spaces that can be used in real-world applications, 
e.g. besides of examples described in Section 3, general formulas of Theorem 2.1 can also be applied in the 
exponential space ^jl, e Al “, e As! “,..., e A "“ [a,0 < Ai < A 2 < ... < A n , a < /3, or in the space 

(^jl, u Xl , u A2 ,..., u Xn :uG [a,/?]|^, 0 < Ai < A 2 < ... < X n , [a,/3] C (0, 00 ) of restricted Muntz polynomials 
among many others. 

Storing in permanent lookup tables the zeroth and higher order endpoint derivatives of the ordinary EC basis 
(1) and of the normalized B-basis (4) induced by it, general formulas (10)—(11) and the proposed control point 
based curve/surface modeling tools can efficiently be implemented up to n = 15, even by means of a sequential 
algorithm. If one uses multi-threading, the value of n, for which one can provide an efficient implementation, 
can be higher. For arbitrarily large values of n, the presented results are mainly of theoretical interest. 
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A Mapping Bernstein polynomials to the ordinary power basis: revisited 

In what follows, we show through direct calculations that the substitution of derivatives (45)-(46) into formulas 
(10)— ( 11) leads to the classical transformation matrix (47) that maps the Bernstein polynomials of degree n to 
the ordinary monomials. In order to perform these direct calculations, we will use the following two lemmas. 
(Naturally, by using simple algebraic manipulations, one can obtain this basis transformation in a much easier 
way. However, we thought it would be interesting to revisit this property in the general context of Theorem 
2 . 1 .) 


Lemma A.l (Sum of multinomial coefficients over proper partitions; [14, Theorem 2.2]) For all p > q 

we have that 


E 


p'- 

iiLi**' 




(55) 


where s q , p is a proper q-fold partition of p, i.e., s q , p is an ordered set of non-negative integers {k z } q z=1 such that 
"Yl q z= i k z = p and kp > 1, W = 1, 2,..., q. (Sq tP denotes the set of all proper q-fold partitions of p.) 


Lemma A.2 For all orders p > 1, one has that 
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From hereon, a number in paranthesis above the equality sign indicates that we apply the corresponding 
identity. For example, when 2 < j < [ and i = 1, 2,... ,j — 1, one has that 
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The remaining cases can be handled in a similar way. 


B A simple Matlab example 

In order to ease the reviewing process, Listing B.2 provides implementation details, by means of which one can 
verify the numerical values listed in Example 3.5. Note that, by using effective multi-threaded object oriented 
C++ programming and OpenGL rendering techniques, even general cases can similarly be treated. With this 
illustrative Matlab code we opted for simplicity. Listing B.2 is based also on Doolittle’s LU decomposition 
which is implemented in Listing B.l for the sake of convenience. 


Listing B.l: LU factorization of a regular real square matrix by means of Doolittle’s algorithm 


7. it is assumed that the real square matrix W is regular 
function [L , U] = Doolittle(W) 

[m, n] = size (W) ; 
if (m ~ = n) 

error (’Matrix must be square!’); 

end 

U = zeros (m); U(l, :) = W(l, :); L = eye (m); 

for i = 2:m 

for j = 1:m 

for k = 1:i-1 
s = 0; 
if (k == 1) 
s = 0; 
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else 

for 

end 

end 

L (i, k) 

end 

for k = i:m 
s = 0 ; 
for p = 
s = 

end 

U(i, k) 

end 


p = 1 : k -1 

s = s + L(i, p) * U(p, 
= (W (i, k) - s) / U (k, 

1 : i -1 

s + L(i, p) * U(p, k) ; 
= W(i, k) - s; 


k) ; 
k) ; 


end 

end 


Listing B.2: Transforming the exponential trigonometric normalized B-basis ( ) to the ordinary basis ( ) 


/ the dimension is n+ 1 = 5, while the user defined growth rate equals uj = -^ 
clear; dim = 5; omega = 1/(3 * pi); 

l fix a proper shape parameter /3 > 0 by using [ 6 ], i.e., ensure that the underlying vector space 
/ spanned by the ordinary basis {<^ 4,1 (u) :uE [0, /3 ]}^_ 0 of the form (49) is indeed an EC space 
beta = 5 * pi / 6; 

l auxiliar symbolic function definitions for symbolic differentiations 
syms u; psi = cell(l, dim); 

psi{2Xu) = exp(-omega * u) * cos(u); psi{3}(u) = exp(-omega * u) * sin(u); 
psi{4Xu) = exp ( omega * u) * cos(u); psi{5Xu) = exp ( omega * u) * sin(u); 

°/ 0 function handles of derivatives (u) : u E [0,/3]j' ; note that array indices start from 1 

phi = cell (dim ); 

for j = 0 : dim -1 

% derivatives of the constant function ipo(u) = 1 have to be handled carefully 
phi-Cj+1, 1 } = @(u)ones(l, length(u)) * (j == 0 ); 

% speed up the calculations by converting all auxiliar symbolic expressions to function handles 
for i = 2 : dim 

phi-Cj+1, i} = matlabFunction(diff (psifi}, u, j)); 

end 

end 


/ define function handles for particular integrals {v± : i (u) : u E [0,/3]}^_ o of the form (38) 
v = cell (dim ); 

for j = 1 : dim 

for i = 1 : dim 

v{j , i} = @(u, rho)rho(i, :) * [phi-Cj, lXu) ; phi{j, 2 }(u); phi-Cj, 3}(u); ... 

phi-Cj , 4}(u); phi{j , 5}(u)]; 

end 

end 

°/ 0 calculate the transformation matrix [pi,j] i L 0 J=0 that maps the ordinary basis (u) : u E [0,/3 ]}^_ 0 

/ to the bicanonical particular integrals {^4 ,i (u) : u E [0,/3 ]}^_ 0 that fulfill the initial conditions (39) 
rho = zeros(dim); 

/ different initial conditions are associated with different extended collociation matrices M 
M = zeros (dim ); 

for i = 1 : dim 

/ current initial condition 
c = zeros (dim , 1 ); c(i) = 1 ; 

7, extended collocation matrix 
for j = 1 : i 

M (j , : ) = [phi{j, 1 X 0 ) , phi-Cj, 2 X 0 ), phi{j, 3X0), phi{j, 4X0), phi-Cj, 5X0)]; 

end 

for j = 1 : dim -i 

M(i + j, :) = [phi{j, lXbeta), phi-Cj, 2}(beta), phi-Cj, 3}(beta), ... 

phi-Cj, 4}(beta) , phi-Cj, 5Xbeta)]; 

end 

l coefficients corresponding to the current initial condition 
rho (i , :) = M \ c ; 

end 

*/. calculate the Wronskian matrix 4 v 4 3 v 4 2 ^4 1 v 4 0 ] (^) reverse£ i system {^ 4,4 — i ( u ) : u E [0,/3 ]}^_ 0 at u = (3 

W_reversed_beta = zeros(dim); 

for i = 1 : dim 

for j = 1 : dim 

W_r e ver sed_bet a (j , i) = v{j , dim + 1 - iXbeta, rho); 

end 


end 
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l perform the Doolittle LU-i actorization of WV i (8) 

r [ v 4,4 > w 4,3 < v 4,2 > v 4,l ’ X7 4,0J V ' 

[L, U] = Doolittle ( W_reversed_beta) ; 

/ calculate the essential parts of the inverse matrices U~ 1 and L -1 
mu = U \ eye(dim, ceil(dim / 2)); 

lambda = L \ eye (dim, 1); 

l define function handles for the zeroth and higher order derivatives of the normalized 
°/ 0 B-basis { 64 , 2 ( 14 ) : u E [0,/3]}^_ 0 ; do not forget the reflection invariant property of the 
°/ 0 underlying EC space 
b = cell (dim ); 


j = 

1: dim 













b{j , 

, 5} = 

0(u) lambda (1, 

1 ) 

* 

mu (1 , 

1) * v{j, 

6}(u, 

rho ) ; 






b{j , 

1} = 

<§(u)(-l)~(j - 

1 ) 

* 

b{j , 

5}(beta - 

u) ; 







b-C j , 

, 4} = 

0(u) lambda (2, 

1 ) 

* 

(mu (1 

, 2) * v{j 

. 5}(u j 

, rho) 

+ mu (2 , 

2) 

* , 

4} (u, 

rho ) ) 

b-C j , 

, 2} = 

®(u)(-l)“(j - 

1 ) 

* 

b{j , 

4}(beta - 

u) ; 







b{j , 

, 3} = 

0(u) lambda (3, 

1 ) 

* 

(mu (1 

, 3) * v{j 

. 5}(u, 

, rho) 

+ mu (2 , 

3) 

* , 

4} (u j 

rho ) 


+ mu(3, 3) * v{j , 3}(u, rho )); 


end 

°/ 0 from hereon, one can use the normalized B-basis functions { 64,2 (u) : u E [0,/3]}^ =0 , 
°/ 0 e.g. we can evaluate and plot them 

figure; hold on; axis equal; color = hsv(dim); u = linspace( 0 , beta); 

title (’Exponential trigonometric normalized B-basis functions of order 4’); 

for i = 1 : dim 

plot(u, b{l , i}(u), ’Color’, color(dim + 1 - i,:)); 

end 


l or we can also create function handles for the derivatives of EC B-curves of the type (17) 
B_c = @(u, p, j) p * Cb-Cj , l}(u); b{j , 2}(u); b{j , 3>(u); b{j , 4>(u); b{j , 5>(u)]; 

/ where p denotes a real matrix of the size <5 X (n + 1) representing the control points 

°/ 0 using Thorem 2.1, one can calculate the entries of the transformation matrix that maps the 
% normalized B-basis {£> 4,2 ( u ) : u E [0,/3 ]}^_ 0 to the ordinary one {<£> 4,2 ( u ) : u E [0,/3 ]}^_ 0 
T = zeros (dim ); 


T (1 

, :) 

= ones 

(1, dim ) ; 



for 

i = 

2 : dim 





T ( i , 

1 ) 

= phi{l, 

i>(0); 



T ( i , 

, 2) 

= phi {1 , 

i}(0) + phi {2, 

i}(0) / b{2 , 2X0) ; 


T ( i , 

, 3) 

= phi {1 , 

i}(0) - phi {2, 

i>(0) / b{3, 3X0) / b{2, 2X0) * b{3, 2X0) 




+ phi{3 

, i>(0) / b{3, 

3X0) ; 


T (i , 

, 4) 

= phi {1 , 

i}(beta) + phi{2, i}(beta) / b{2 , 4}(beta); 


T ( i , 

dim) 

= phi{1, 

i} (beta) ; 



end 


l after this, one can apply Theorems 2.2—2.3, or Algorithm 2.1 and its extension to obtain 
/ possible control point configurations; e.g. assume that one intends to describe the arc 
i cos (u + 7 ) , e u ^ u '*’ 7 ^ sin (u + 7 )] T , w E [0, of a logarithmic spiral, where 7 EK is an 

l arbitrarily fixed phase change 

figure; hold on; axis equal; u = linspace(0, beta, 30); 

title (’Control point based exact description of different arcs of a logarithmic spiral’); 
gamma = (-dimrdim) * beta; m = length ( gamma ); color = hsv(m); 

for i = 1 : m 

/ calculate the position of control points corresponding to the ith phase change 

p = [exp(omega * gamma(i)) * (T(4,:) * cos ( gamma (i)) - T(5,:) * sin ( gamma (i))) ; ... 

exp(omega * gamma(i)) * (T(4,:) * sin ( gamma (i ) ) + T(5,:) * cos ( gamma ( i )))] ; 

l evaluate and render the points and tangents of the EC B-curve generated by p 

arc = B_c(u, p, 1); hodograph = B_c(u, p, 2); tangent = arc + hodograph ; 
plot ( arc (1,:), arc (2,:), ’Color’, color(m + 1 - i,:), ’LineWidth’, 4); 

plot([arc(l,:); tangent (1,:)] , [arc (2 , :); tangent (2,:)], ’Color’, color(m + 1 - i , : )); 
# /. render the control polygon 

plot(p(l,:), p(2,:), ’-o’, ’Color’, color(m + 1 - i,:), ’LineWidth’, 2); 

end 


















